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We address the importance of symmetry and symmetry breaking on linear response theories 
of fermionic BCS superfluids. The linear theory of a noninteracting Fermi gas is reviewed and 
several consistency constraints are verified. The challenge to formulate linear response theories of 
BCS superfluids consistent with density and spin conservation laws comes from the presence of a 
broken U(1)em symmetry associated with electromagnetism (EM) and we discuss two routes for 
circumventing this. The first route follows Nambu's integral-equation approach for the EM vertex 
function, but this method is not specific for BCS superfluids. We focus on the second route based on 
a consistent-fluctuation-of-the order-parameter (CFOP) approach where the gauge transformation 
and the fluctuations of the order parameter are treated on equal footing. The CFOP approach 
allows one to explicitly verify several important constraints: The EM vertex satisfies not only a 
Ward identity which guarantees charge conservation but also a Q-limit Ward identity associated 
with the compressibility sum rule. In contrast, the spin degrees of freedom associated with another 
U(l)z symmetry are not affected by the Cooper-pair condensation that breaks only the U(1)em 
symmetry. As a consequence the collective modes from the fluctuations of the order parameter only 
couple to the density response function but decouple from the spin response function, which reflects 
the different fates of the two U(l) symmetries in the superfluid phase. Our formulation lays the 
ground work for application to more general theories of BCS-Bose Einstein Condensation (BEC) 
crossover both above and below Tc. 

PACS numbers: 74.20.Fg,74.25.N-,03.75.Ss 

I. INTRODUCTION 

Linear response theories have been an important tool for studying transport and dynamic properties of many-body 
systems [TH3]. Although there have been myriad successful examples in classical or non-interacting systems, developing 
a linear response theory for complex systems or strongly correlated systems could be quite a challenge. As summarized 
in Ref. [3], it is difficult to obtain consistent expressions of the compressibility of an interacting electron gas from 
the derivatives of thermodynamic quantities and from the correlation functions. A naive calculation of the response 
functions of BCS theory of conventional superconductors could have led to a violation of the charge conservation. 
The spontaneous symmetry breaking of the superfluid phase further complicates the treatment of any linear response 
theory of BCS superfluids. Thus it requires more sophisticated treatments [Hl^ to obtain the response functions that 
respect conservation laws. 

Here we analyze two consistent linear response theories of BCS superfluids in great details and review the importance 
of gauge invariance, the (generalized) Ward identity, sum rules, and how they impose constraints on response functions. 
Also crucial is an additional Q-limit Ward identity [6] which is related to the compressibility sum rule. The non- 
interacting Fermi gas satisfies all constraints, albeit in a trivial sense. In the presence of interactions, we will focus 
on a linear response theory which we call consistent-fluctuation-of-the order-parameter (CFOP) theory for fermionic 
superfluids. Importantly, at the mean-field level, this approach is consistent with all conservation laws and associated 
sum rules. This will set up a solid foundation for the generalization to the more general theories of BCS-BEC 
crossover 0. We also present a spin linear response theory in a systematic manner to demonstrate that fundamental 
constraints including the sum rules and Ward identities are also satisfied in the spin channel. One major difference 
between the density and spin response functions is that in the superfluid phase the order parameter, which corresponds 
to the condensate of Cooper pairs, breaks only the symmetry associated with the density response functions, but not 
the symmetry associated with the spin response function. Therefore the density response function exhibits richer 
structures across the superfluid phase transition. This is the base for constructing an order-parameter like quantity 
that measures the "spin-charge separation" in a BCS superfluid as proposed in Ref. [8]. 

There have been two major approaches for addressing the linear response of fermionic BCS superfluids. The first one 
is attributed to Nambu [4| , who reformulated the problem in the two-dimensional Nambu space and pointed out that 
if certain types of corrections to the EM interacting vertex are considered consistently with those corrections to the 
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self-energy in the Green's function, the gauge invariance can be maintained exphcitly. Nambu presented a generalized 
form of Ward identity (GWI) and proposed an integral equation for the full electromagnetic (EM) interacting vertex 
in the Nambu space. He stated that this EM vertex must satisfy the GWI without giving a proof and we present 
our own proof in Appendix [C] By solving this integral equation at small frequency and momentum limits, it was 
shown that the excitations of collective modes correspond to the poles of the density response function. Hence the 
many-particle effects are indeed included in the corrections of EM vertex. Despite many virtues, this method is 
relatively difficult to implement. Firstly, it is very hard to solve the integral equation for the vertex. If one truncates 
the integral equation, one may not get a gauge invariant solution and can not even find the correct collective-mode 
excitations. Secondly, this approach was formulated in the Nambu space and it is difficult to translate the results 
in the two-dimensional Nambu space to their counterparts in the one-dimensional representation space of fermion 
operators. This limits its applicability to more general problems such as BCS-BEC crossover. Moreover, it can be 
shown that the GWI discussed in Nambu's original paper |4| is not a unique constraint for gauge invariance of BCS 
theory. In other words, Nambu's method is not specific to BCS theory: It is more general but harder to obtain an 
exact solution. 

There is yet another approach which we call the CFOP approach with a totally different structure when compared 
to Nambu's method. Our goal is to explain this approach and test its results against some fundamental constraints. 
Within this approach, we treat the gauge transformation and the fluctuations of order parameter on equal footing, 
such that many-particle effects are also explicitly included in the EM vertex. This approach was first proposed by 
Kadanoff and Martin[5, in a less complete form. They only considered the fluctuations of the order parameter, and 
tried to decompose the three-particle Green's functions in a way that can respect gauge invariance. It was then 
independently formulated in several unrelated papers by Betbeder Matibet and Nozieres [9j and Kulik et al. [10] 
in a more complete form where both phase and amplitude fluctuations are considered. Later on, Zha and Levin[TT] 
revisited it and presented three identities of response functions, which are now known to be part of the Wis for 
response functions. Recently, it was again proposed in the Keldysh formalism with time-ordered Green's functions 
in Ref. [12]. A similar derivation using a kinetic-equation formalism is also discussed in Ref. [13]. In this paper, we 
cast this formalism in a more systematic and covariant form. We will present more virtues of this approach, such as 
the consistency with the Q-limit Ward identity and the /-sum rule. The EM vertex and its generalized form can be 
both found in different representation spaces. A comparison between Nambu's method and the CFOP theory will 
also be presented. We will show that the CFOP approach is indeed a consistent and manageable linear response 
theory for BCS superfluids. The generalization of the CFOP theory to relativistic BCS superfluids can be found in 
Ref. jl4| . Besides linear response theories discussed here, one may flnd more discussions on the gauge invariance of 
BCS superconductors in terms of effective field theories in Ref. |15l . 

This paper is organized as follows. In Section [H] we briefiy discuss the general symmetries in the theory of two- 
component Fermi gases with contact interactions. In Section |III| we review the density and spin linear response 
theories for non-interacting Fermi gases. In SectionjlV] we spend a significant part of the paper on the CFOP linear 
response theory of BCS superfluids. More speciflcally. Section |IV A| presents the derivation and some results of the 
CFOP linear response theory; Section IV B addresses the (generalized) Ward identity and Q-limit (generalized) Ward 
identity of the CFOP linear response theory; Sections |IV C| and |IV D| review Nambu's integral-equation approach and 
present a comparison of our CFOP approach with Nambu's method. In Section [V] we develop a parallel formalism for 
the linear response theory in the spin channel by "gauging" the U(l)2 symmetry. Section VI concludes our work. 



II. SYMMETRIES OF THEORIES OF FERMI GASES WITH CONTACT INTERACTIONS 



Throughout this paper, we follow the convention e = c = h = 1 and use a to denote the spin or pseudo-spin t , i 
with t and a being the opposite of J, and a respectively. The metric tensor of the Minkowski spacetime is chosen as 
rj^"^ = diag(l, —1, —1, —1) and Einstein summation convention is adopted. For a two component Fermi gas interacting 
via contact interactions, we consider the Hamiltonian 

J d'^^l{^)(^^^t?jM^)~g J d3xV;|(x)^|(x)V'i(x)V't(x), (1) 

where tp and tp^ are the annihilation and creation operators of fermions, n is the chemical potential, m is the fermion 
mass, and g is the bare coupling constant. There is an implicit summation over the pseudo-spin indices a. The 
Hamiltonian ([T]) has a SU(2)xSU(2) symmetry [HlinillT]. The flrst SU(2) symmetry is generated by 

— >■ coshxV't + sinhxV'l : i'l sinhx^J -|- coshx"!/"-!- , 
— coshx?/;^ — isinhxip] i "0! ~^ isinhxi/ij 4- coshx^/;^ , 
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^ e-'^^„, ^ e'^^l (2) 

where x is a continuous parameter. The transformation on the third hne is the weU-known U(l) symmetry associated 
to electromagnetism if the system is charged. The generators of the these transformations are —iai, —ia2 and (73 
in the space spanned by the Nambu spinor {ip^ , tp'^^)'^ , where CT; are the Pauh matrices. Hence the symmetry ^ is 
indeed SU(2), or more precisely, SU(1,1) with the U(1)em being its subgroup. 
The second SU(2) symmetry is given by 

tp-f — >■ cosa'i/'t + jsina'0^, ^ —isinatpJ^ + cosai/'J, 
'0^ — cosa-0^ + s'maipi, i/j^ — > — sinaV'l + cosaV"!, 

e-''^-"V<T, V'i ^ e'^^-^^jl, where 5^,^ - ±1. (3) 

The generators are cti, 172 and CT3 in the space spanned by {•ip^jipi)'^ . The transformation on the third hne is the spin 
rotation around the z— axis, which forms the subgroup U(l)z of the second SU(2). In what follows, we will focus on 
the U(1)emxU(1)z symmetries, and we may call the theories associated with them as being in the density and spin 
channels respectively. 

When the continuous parameter x become space-time dependent, the two global U(l) symmetries are "gauged" 
respectively. To keep the Lagrangian invariant under the gauge transformations, one needs to couple the fermion field 
by an external vector field A^^ — (0, A) which transforms as A A — Vx, 'Z' ^ + For the U(1)em symmetry, 
the Hamiltonian (24) becomes 



:V4(x)( ^P 2m / d'M^)^U^)M^) - 9 j rf^xV4(x)V'|(x)0^(x)0^(x) 

J d^x^t (x)(|^ „ M)V'a(x) + J d3xJ(x) • A(x) + J d3xn(x)0(x) -gj d^x^j (^)^t (^)^^(^)^^(^)^ (4) 



where 



J(x) 



1 

2 mi 



V4(x)(V^,(x)) - (VV4(x))^.(x) - -A(x)V'i(x)0,(x), 



1 



m 



n(x)=0t(x)v,„(x). (5) 

Here n(x) is the particle number density and J(x) is the mass current. There is also an implicit summation over 
the repeated Greek index a. When the external field is the EM field, the Noether current for the global U(1)em 
symmetry is = in,J) which obeys the conservation law 9^J^ = in the Heisenberg picture, where 9^ — (^, V). 

The same discussion can be implemented in the spin channel. The fermion can also couple to an effective external 
field A^. Hence the Hamiltonian with a U(l)z "gauge" symmetry is given by 

I d^x^'-T'"^' - ^-A(x)) 



V;i(x)(^^^^-A^-Ai)V,(x)+ J d^^<f>{^)S,^U^)ijA^) - g J d3xV4(x)V'|(x)V^(x)V't(x) 
d^x^At (x)(|^ - /i)V'.(x) + J d3xJs(x) • A(x) + J d3xns(x)0(x) ~gj d3xV;|(x)V'|(x)V'i(x)7A^(x), (6) 



where 



Js(x) = -^5,[v4(x)(VV'.(x)) - (V0t(x))^^(x)] _ lA(x)Vi(x)0^(x), 

ns(x) = 5,Vi(x)V'.(x). (7) 

Here S^^i — ±1 as in The external field has a different physical meaning from that in the density channel. 
Since ng is the z component of spin density, the field (p coupled to ng corresponds to . Js is the difference between 
the spin-up and the spin-down currents, i.e., the magnetization current. Therefore the field which couples to it is 
the magnetizing field A = m. The effective external vector field is thus = (i?^ , m) . The resulting Hamiltonian 
then describes a generalized spin-magnetic field interaction. The Noether current for the global \J{l)z symmetry is 
Jg = {ns, Js), which also satisfies the conservation law 9;^ Jg = 0. 
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III. LINEAR RESPONSE THEORY FOR NON-INTERACTING FERMI GASES 

We begin with the non-interacting Fermi gases where g ~ 0. Many important identities can be verified exphcitly 
and they provide useful hints for the development of consistent theories of linear response for interacting Fermi 
gases. When a non-interacting two-component Fermi gas is perturbed by a weak external gauge field, one can discuss 
the response in the density channel as well as in the spin channel the leading-order or linear-approximation of the 
perturbation. 

In momentum space, the linear response theories associated with the two U(l) symmetries mentioned above can be 
expressed in a unified form as H = Hq + g , where 

P 

Hd=Y1 ^ip+q7^(P + q> P)A^,q^ap, = X! ^ip+q7s^(P + '^^ P)^MqV'<Tp- (8) 

pq pq 

Here Cp = 2m ~ ^^"-^ have introduced the interacting vertices 7'^(p + q, p) — (1, for the density channel 

and 7go.(p + q, p) = {Sa, S^^^^^^) for the spin channel. The latter explicitly depends on the spin (or pseudo-spin) 
and Figjl] illustrates this spin dependence. The four-currents of the density and spin channels in momentum space 
are given by 

"^q = E ^ip^''(P + P)^-P+q' -^Sq = E ^U^sAP + P)^.p+q- (9) 
P P 

Note in the perturbed Hamiltonian g , we only keep the terms up to the linear order of the external field and 
higher order terms are neglected. Hence in the expressions of the currents and Jg^ the linear terms of are 
dropped, which is different from the currents given by Eqs. ([5]) and Q. 

For the rest of this section we will focus on the density channel since the linear response theory in the spin channel 
has a similar structure as the former, and the response kernels can be shown to be the same. The imaginary time 
T = it is introduced here and the Heisenberg operator is defined as 0{t) = e^'^Oe~^'^ . The spacetime coordinates 
are x = (r, x). Hence the Green's function is defined by Go{x, x') — — (T^[^CT(a;)V'J(^')]) where denotes the r-order 
of operators. Its expression in momentum space is Go(«w„,p) = (iwn — ^p)~^. The four-momentum is defined as 
P = = p) where ojn = {2n + \)nkBT is the fermion Matsubara frequency. The particle density is given by 

n = r^^Go(P) = 2^/(ep), (10) 

iU)n P.CT P 

where f{x) = \ j[e^l^'^'^ + 1) is the fermion distribution function. By defining h^^ = —ri'^'^{l — rj'^^), we formally write 
the perturbation of the EM current J'^ as 

6J^{t, q) = 5](^ip(T)7'^(p + q, p)V<xp+q(r)) + -h^^'' A,{t, q), (11) 
p 

The extra term ■^h^'^ of 5J^^ arises from the second term of J or Jg (See Eqs.([5]) or ([t])), since J or Jg already contains 
a first order term of A^ . The linear response theory is then written in the form 

<5J^(r,q)- j dr'{Q^^%T-T',ci) + ^h^^5{T-T'))A,{T\c^), (12) 

where the response kernels are 

Qr(T-r',q) = -(T,[J^(r,q)J''(r',-q)]) 

= -Y^iTMpirh'-iv + q, p)^.p+q(T)Vl,p,+q(T')7^(p', p' + q)^.'p'(T')]). (13) 
pp' 

Implementing a Fourier transform and making use of Wick's theorem, we obtain 

Q^"^(*r!z, q) - 2r^ ^ 7^(p + q, p)Gop+q,p'+q(zc^„ + 1^1)^ {p\ p' + q)G'op,p' (zw„) 

iuJn pp' 
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E ^"(^ + ^)Go(P + Qh^iP, P + Q)Go(P). (14) 



where Q = — (if2;,q), fi; is the boson Matsubara frequency, and Gop,p'(jwn) — Gop(i^n)'5p,p' with Gop(cLi„) = 
Go{P). For convenience, we have defined 7''(P + Q, P) = 7^(p + q, p)- The factor 2 comes from the summation 
over the spin (or pseudo-spin) indices a. The spin response kernels have the same structure as one can see from the 
facts that ^g„{p + q, p) = S'cr7''(p + q, p) and = 1. Note there are also two spin interacting vertices inside the 
expression of the spin response kernels similar to (13 1. The detailed expressions of the response functions are listed 
in Appendix |Aj where = £,p± a . 

Due to the simple structure, one may verify the following identities explicitly: (1) Ward identities (Wis), (2) 
compressibility sum rule, (3) /-sum rule, and (4) Q-limit Ward identity. By direct evaluating each term, one can 
verify that the EM and spin interacting vertices satisfy the Wis 

q^nP + Q,P) = G^\P + Q)-G^\P), 

Q^lLiP + Q^P) = S,{G^\P + Q)^G^\P)). (15) 
This leads to the Wis for the response kernels q^Q'^'^ [Q) — where Qq'' — Q'^'^ + m^^'^- shown as follows 

(Q) = 2 V[Go i(P + Q)- Go\P)]GoiP + Q)Y{P,P + Q)GoiP) - -g^il - v"") 



p 

2Y,Go{PW{P + Q,P)- Y{P - Q, P)] - -q'ii - v"") 



m 



2^] Go(P)^(i - v""') - -q^a - v''') = 0, (16) 



where = PX^ii^ Wis for the response kernels further lead to the conservation of the perturbed current 

q^5J^{Q) = 0. For the spin response, the same conclusions can be obtained. 

Next we show that the response function satisfies the compressibility sum rule [T and /-sum rule 

(til 

— = -g™(c.-0,q^O), (17) 

dwwxpp(w,q) = , (18) 

where Xpp — ~:^I™'9o° is the density susceptibility and the analytical continuation i^i — )• a; + iQ^ has been applied. 
Although these two sum rules can be directly proven by the explicit expressions of the response functions given in 
Appendix |Aj here we give a more instructive proof that has a nice connection with the U(l) symmetries in the 
Hamiltonian. For the compressibility sum rule, we have 

|l = 2 E = -2 E Go(^)^^|P = -2 Imr J] Go(P + 0)Go(P)U.o = -0S°(- = 0, q ^ 0), (19) 



where the expression ( 14) and 7'^(P + Q, P) = 1 have been applied. In fact, this is a special case of the Q-limit Wl 
r° = 1 — {d'S/d^), where F" and 7° are the full and bare vertex functions, S is the self energy of fermions, and the 
limit w = 0, q — has been taken. For non-interacting Fermi gases S = and r° = 7° = 1 so the Q-limit WI is 
trivially satisfied. For the /-sum rule, we need to implement the real time formalism to describe the non-equilibrium 
transport theory. The real time response function corresponding to Eq.(13) is 

- t': q) - -^Ht - t'){[J^{t,q), r{t\ -q)]), (20) 

where the time-dependent Heisenberg operator is defined by 0{t) — e^^'Oe^'-^*. The imaginary part of Q^'^ is 

imQ^ - 1', q) = -^([^^(i, q), •rit'. -q)])- (21) 

Using ujIthQq'^ = 1u\{ujQq'^)^ we have 

1 f^°° dn(t q) 
/ dujujlmQ'^'"{uj,q) ^ {[i —] — , n(i, -q))^ = ([q • J(t, q), n(i, -q))<^ 
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Figure 1: The diagrams for the spin vertex function. Different spin indices have difference signs. 



= • ^^(e'""* WpW^P+q^i'p'+qW.'p' - V'I.p.+q{^ip,^.'p'}^.p+q)e-^^*). 



PP 



(P+ 2) - (P- 2)/ ,t n\ I M\ 

p 



(22) 



where in the first hne we used the conservation of the current 9^J'^ — 0, in the last hne we changed the variables 
p — >■ p — q for the second term. The subscript uj means the Fourier transform of the function with the argument tu. 
One thus see that the two sum rules are closely connected to the U(1)em symmetry. 



IV. GAUGE-INVARIANT LINEAR RESPONSE THEORY FROM THE CFOP APPROACH FOR 

FERMIONIC SUPERFLUIDS 



In BCS theory of fermionic superfluids, the order parameter is given by 

A(x) =g(V.^(x)7/^^(x)>. (23) 
Following the mean-field BCS approximation, the Hamiltonian without any external field becomes 

H = J d3x^t(^)(^|!_ _ ^)^^(x) - J d3x(A(x)^:f(x)^J(x) +h.c.). (24) 

One may see that the U(1)em symmetry is spontaneously broken if A(x) 7^ below Tj, while the U{l)z symmetry 
remains intact. Hence the original Hamiltonian ([T]) has a U(l)xU(l) symmetry while the BCS Hamiltonian ([24]) in 
the broken-symmetry phase only has a U(l) symmetry. The phase with a broken U(1)em symmetry below Tc brings 
challenges of how to cast its associated linear response theory in a gauge invariant form. Below Tc, the breaking of 
U(1)em and the unbroken U(l)z symmetry may be called a "spin-charge separation" in BCS theory. We will explore 
this phenomenon in depth after we present the consistent linear response theories in both density and spin channels. 
Above Tc, both symmetries are not broken and the charge and spin degrees of freedom are not separated there as 
they do below Tc. 

Here we present a U(1)em gauge-invariant linear response theory for BCS superfluids based on the consistent 
fluctuation of the order-parameter (CFOP) approach. The gauge invariance is basically a one-particle problem from 
the point of view of quantum field theory while BCS theory is essentially a many-body theory (for a review, see jl8l ). 
This contrast highlights the importance of Ward identity (or its generalized form) since one can check the gauge 
invariance of a theory by verifying the corresponding WI and we will give some concrete examples. 



A. Nambu Based Notation for Linear Response 



As we have seen previously, there are terminologies like the Nambu space, one-dimensional space and the "Ward 
identity" (WI) associated with different physical quantities. Here we explain them in details. The Nambu space is the 
space in which Nambu spinors are defined. For non-relativistic theories, it is a two-dimensional space and operators 
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are written as two by two matrices. The basic framework of the CFOP hnear response theory for BCS superfluids 
is easier to develop in the Nambu space. The terminology one-dimensional space is the abbreviation of the one- 
dimensional representation space of fermion operators. The representation space of fermion operators is actually the 
four-dimensional spinor representation space of the Lorentz group. In the non-relativistic limit one may approximate 
it by a one-dimensional space representation if the spin is simply labeled as a subscript. Most linear response theories 
for normal Fermi gases are formulated in this space since the corrections to the interacting vertex between fcrmions 
and external fields are hard to be presented in the Nambu space. The Ward identity in quantum field theories refers to 
the relation among the response (correlation) functions that describe the effects of symmetry transformations allowed 
by the theory [T^ . Here we use the Ward identity specifically for a diagrammatic identity between the vertex function 
and the fermion propagator, which reflects the symmetry from the EM gauge transformation. In the Nambu space, we 
will adopt Nambu's convention to characterize each of these identities as a "generalized Ward identity" (GWI). There 
are seminal reviews [H [TH] on how to cast BCS theory in the Nambu space. However, there are less reviews on the 
same approach in the one-dimensional space [7] . A fully gauge invariant EM interacting vertex in the one-dimensional 
space has not been studied systematically. Later on, we will present one by using the results found in the Nambu 
space. Here we first formulate the linear response theory in the Nambu space for BCS superfluids based on the CFOP 
approach. 

We define <7± — ^{(Ji ± i(T2) in the Nambu space and introduce the Nambu-Gorkov spinors 



p 



4-r> 



vl'^ = [4p>;-p]. (25) 



Here ai, i — 1,2,3 are the Pauli matrices. The Hamiltonian Q in momentum space after the mean-field BCS 
approximation becomes 



H = J2 *pep^3*p + J2 *p+q( - ^^^^ + - Aq^+ - A-q'^-j^p, (26) 



The order parameter ( |23[ ) in momentum space is generalized to include finite-momentum fluctuations from its equilib- 
rium value so Aq = g2!^(*p'^-^p+q)! which will be imposed as a self-consistency condition. When the external EM 
field is weak, the order parameter is perturbed and deviates from its equilibrium value. The order parameter in equilib- 
rium is A, which is at q = and can be chosen to be real by the U(1)em symmetry. We denote the small perturbation 
of the order parameter as A^ so Aq = A A'^. By introducing Aiq = -(A^ + A'l^)/2 and Asq = -^(A;^ - A'l^)/2, 



the Hamiltonian ( 26 1 splits into two parts as H = Hq + H' with one containing the equilibrium quantities and the 



other containing the deviation from the equilibrium. 

i/o = 5I*p^P*P' ^'-E*p+q(^iq^i + ^2qa2-l-A,,q7'^(p + q,p))*p, (27) 
p pq 

where E-p = ^pcra — Acti is an energy operator and 7^(p -I- q, p) — (cs, is the bare EM vertex in the Nambu 

space. As discussed in Ref. [10], the fluctuations Aiq and A2q introduce the amplitude mode and phase mode as the 
collective modes introduced by the perturbation of the order parameter. 

For the equilibrium par, the quasi-particle energy is given by Ep — + IS?. The propagator in the Nambu space 

is 

<S(P).6.(i.„)^-^^(?<JJ _1\%X (28) 



where 



G{P) = -. V+ ■ ' ^(^) = -"p^p(. ^- . ) (29) 

lUJn - Ep lUJn + Ep \lUJn - Ep lUJn + Ep 



are the single-particle Green's function and anomalous Green's function respectively and Wp,fp — ^(1 ± If")- 
The number density and the order parameter in equilibrium can be expressed by the propagator as 



n = Tr^(a3G(P)), A = gTr ^ (aiG(F)) , 
p p 



(30) 
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which give the number and gap equations 

We can also introduce the counterpart of the coefficients u,v in the Nambu space as Up,Vp — ^(lzb|^). Using the 
properties Up — Up, = Vp and UpVp = 0, one can show that in the Nambu space the propagator has the form 

U V 

G{P) = , L + . L. (32) 

liOn - Ep ILOn + Ep 



This has a similar expression as the single particle Green's function shown in Eq. ( |29[ ). In fact, any function of 
iujn — Ep, for example F{iLUn — Ep), can be expressed as F(ia;„ — Ep) = Up{iujn — Ep) + Vp{iuj„ + Ep) [T3]. It is an 
interesting property of BCS theory in the Nambu space and will be useful for deriving the response functions. 
The interacting Hamiltonian H' in Eqs. ([27| can be cast in the form 



^' = E*P+q*q-^(P + q'P)*P' (33) 

pq 

where 

*q= (Alq,A2q,A^q)^, £ (p + q, p) = (cti , CJa , 7^ (p + Q, p) ) ^ , (34) 

are defined as the generalized driving potential and generalized interacting vertex, respectively. By using the imaginary 
time formalism, we assume that under the perturbation from H' , the generalized perturbed current SJ is given by 



6J{r, q) = V(*t,(T)S(p + q, p)*p+q(T)) + A,{r, q). (35) 



Here JJg denotes the perturbed EM current and 6Ji^2 denote the perturbations of the gap function. The linear 
response theory is written in a matrix form 

(5J(w,q) = 0(a;,q) • *(w,q) 

Qii(w,q) Qi2(a;,q) Q\^{uj,ci) \ / Ai(a;,q) 

g2i(c^,q) Q22{io,q) Qhi^,ci) A2(w,q) |. (36) 

g^\(c.,q) Q^'2(^,q) Qg(w,q) + ^/i^V \A,ii,,q) 

The response functions Q.^ are 

0,,(r - r',q) = - J2{Tr[¥^{r)±^,{p + q, p)*p+q(r)vl/^,^q(T')S,(p', p' + q)*p'(^')]>, (37) 
pp' 

After applying Fourier transform and Wick's theorem, we obtain 

(zf7i, q) = TrT^ ^ (E,(P + Q, P)G{P + Q)t,{P, P + Q)G{P)), (38) 

iUn P 

For convenience, we have defined 7^(P + Q,P) = ^^{p + q, p), hence S(P + Q,P) — S(p + q, p). By using the 
expression ( 32 1 the response functions are evaluated and shown in Appendix [A] 

The perturbation of the order parameter and the EM perturbation are treated on equal footing and this will 
naturally lead to the gauge invariance of the CFOP linear response theory. The gap equation gives the self-consistent 
condition 5Ji,2 — — |Ai,2. Applying this relation to Eq.(36), we get 

^^^_^^3^22^QM}1a,, A2^-^^^^"~^"^^^A.. (39) 

Q11Q22 ~ Q12Q21 '3ii'322 ^ Q12Q21 

where Qii = | + Qii ^nd Q22 = f + ^22- After substituting the results into 

<5J^ = Q^iAi + Q^2A2 + {Q'i,l + -h^^^Au, (40) 
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we get 5J^ = K^^ A^^ where the corrected EM response kernel K^^^ including the effects of fluctuations of the order 
parameter is given by 



Q11Q32Q23 + Q22Q31Q13 ~ Q12Q31Q23 ~ Q21Q32Q13 
Q11Q22 — Q12Q21 



(41) 



Here Q33 = Q33 + ^h'^'^. The gauge invariance condition for the perturbed current, g^fJJ'' = 0, is satisfied if 
q^K^'^{Q) = 0. This is further guaranteed by the following Wis for the response functions 



31 



-2iAQ2i, q^Q 



32 



'2iAQ2 



33 



-2jAQ2"3. 



(42) 



These relations are more complicated than QfiQ^'^ {Q) = for noninteracting Fermi gases since the perturbation of the 
order parameter comes into the linear response theory. The gauge invariance condition of the perturbed EM current 
is immediately derived from Eq. ( 42 ) 



T^^.v o /^mv , r, - f.QllQ'22Q23+ Q'22Q2lQl3- Ql2Q2lQ23^ Q21Q22Q13 



Q11Q22 — Q12Q21 



-2i AQ^3 + 2iAQ^3 = 0. 



(43) 



Now we sketch the proof of the Wis (42). The bare inverse propagator in the Nambu space is Gq^{P) = ii^n^ S.p'^3 
and the self energy is S = — Aui [18]. Hence from Eq.(28) the full inverse propagator can be expressed as G^^{P) = 
Gq"^(P) — S, and one can verify that 



fT3G-i(P + g) - G-\P)a^ - iilin; - (Cp+q - ^p) + 2zA<72 = q^r{P + Q,P) + ^lAa^. 



(44) 



In fact, as we will point out later, this identity is actually the GWI that connects the full EM interacting vertex and 
the full propagator in the Nambu space. Now we show how the Wis are derived from Eq.(44). For the first identity 
of the Wis, we have 



p 

= Tr^ [i<J2GiP)] + Tr^ [G(P + Q)ia2] = 2Tr^ [^'^2G{P)] = 0, 
p p p 

where the cyclic property of the trace has been applied. For the second identity of Wis, we have 

g^Q^2 + 2zAQ22 = Tr^^ [(a3G-i(P + Q) - G-1(P)(T3)G(F + Q)a2G(P) 
p 

A 



(45) 



-2Tr^ [i(7iG{P)] ^ -4iTr^F(P) 
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(46) 



Therefore qfj,Q^2 = — 2iA((522 + |) = —2iAQ22 and we get the second Wl for the response functions. For the last 
WI, we have 

q^Q^3^ + 2zAQ^3 = Tr^[a3r(P + Q,P)G(P)] -Tr^[G(P + Q)r(P,P + Q)^3l - "'/'^(l-'?''") 

p p 

- J2 Tr{G{P)a4r{P + Q,P)- Y{P " Q, P)]) " -q^l " r^") 



TO 



TO 



;i - 77"-") E Tr(G(P)a3) - -q'' {I - v"') = 0, 



(47) 



where the fact that 173 commutes with 7'' has been applied and the number equation ( 30 ) has been used. 
The CFOP linear response theory also satisfies the /-sum rule 



+ C30 



TO 



(48) 
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Here Xpp = — ^ImX"" with K^^ given by the OO-component of Eq. (41 ): 



QllQ%2Q2d. + Q22QziQ% — Q12Q31Q23 ^ Q2lQ32'3l3 



Q11Q22 ~ Q12Q21 

The following lemma will be useful for the proof of the /-sum rule 

+00 2 

da;-Im[q- Q!]3(a;,q)] = n 



(49) 



(50) 



The proof of this lemma can be found in Appendix [b[ From cjlm/v'^*' = Im(ti;A'™) and the third equation of the Wis 
( 42 1 , we have 



Note Q12 



-Q2I, (323 



1 



q • Q°3 - 2i AQ23 

Q11Q22 — Q12Q21 

C (giiq • Q32Q23 - 2zAQnQ22Q§3 + <322q • Q3lQ?3 " 2lAQ22Q2lQ?3 
-Ol2q • Q31Q23 + 2zAQi2Q2lQ23 " Q2iq • Q32Q?3 + 2iAQ2lQ22Q?3)- 



00 



q-Q3c 



Qm^ Q23 = -Q32, (3i3 = Q31 and Q13 = Q31, we have 

Qiiq ■ Q32Q23 + Q22q • Q31Q13 ~ Qi2q ■ Q31Q23 ~ '92iq • Q32Q1, 
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33 



Q11Q22 — Q12Q21 

Using the lemma, it can be shown that proving the /-sum rule is equivalent to proving 



+00 



rfu-'im ' ^-^^^"s + Q22q • Q31Q13 ~ Qi2q • Q31Q23 ~ Q2iq • Q32Q13 _ g 

Q11Q22 — Q12Q21 



(51) 



(52) 



(53) 



Note that Qn, Q221 Q13 and Q23 are even functions of cj, while Q12, Q13 and Q23 are odd functions of w. Hence the 
integrand is an odd function of oj, then the integral indeed vanishes. Therefore, the /-sum rule (48) is satisfied by 
the density linear response theory. From the proof we can see that we need the explicit expressions of the response 
functions. In fact we can not prove the /-sum rule simply by using the Wis (42) or the conservation law of current 
dfj^J'^ = alone, which is different from what we have done in the case of non-interacting Fermi gases. 



B. WI and Q-limit WI of the CFOP Linear Response Theory of BCS Superfluids 

1. GWI and Q-limit GWI in the Nambu space 



To complete our discussions on the CFOP theory in the Nambu space, we now investigate the GWI which connects 
the full EM vertex and the fermion propagator. Due to different representations for the fermions, we will present the 
GWI in the Nambu space as well as the WI in the one-dimensional space One can verify that the bare EM vertex 
and the propagator given in Sec. |IV A| satisfies the bare GWI in the Nambu space 

q^riP + Q,P) = '^sG^'iP + Q) ~ G^\P)'J3. (54) 
Hence, as pointed by Nambu the GWI for full EM vertex has the similar structure 

q^rf^iP + Q,P) = <y3G-\P + Q)- G-\P)a^, (55) 

Moreover, the corrected EM response kernel K'^^ given by Eq. ( [4l] ) should be evaluated by the full EM vertex and 
propagator in the Nambu space as 



X'"'(g) = Tr ^ (f''(P + Q, P)G[P + QW{P. P + Q)G[P)) 
p 



(56) 



a similar form in the one dimensional space will be presented later in Eq.( 73 ). In fact, this full EM vertex Tf^{P+Q, P) 
can be inferred from the expression of isT^". We define 



H^ 



Q31 


Q21 


Q32 


Q22 


Qu 


Q12 


Q21 


Q22 



H^ 



Q32 


Q12 


Q31 


Qii 


Qu 


Q12 


Q21 


Q22 



(57) 
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Hence, from Eq.(41| K^'^ can be expressed as 



Tr ^ {riP + Q, P)G{P + QWiP, P + g)G(P)) - n5'(g)Tr ^ (aiG(F + QW {P, P + Q)G{P)) 



-n^(g)Tr^ {a^GiP + QW{P, P + Q)G{P)) 
p 

Tr^ ([7^(p + g, p) - aini(g) - cT^iv^{Q)]G{p + g)^ + g)G(i^)) 



(58) 



where in the second Une we have substituted the expression (38 1 for Q'^'^, Qi^ and Q23- Comparing Eq.(56) with the 
last hne of Eq. ( 58 1 , one can find that the full EM vertex is given by 



f ^(p + g, p) = riP + Q-P)- oiKiQ) - ^2n'2'(g). 

By applying the Wis ([42|, one can see that I\-i 2 satisfies <7^n^(g) — and q^n2(g) 
vertex given by Eq. (|59[) further satisfies 



(59) 

-2iA. Hence, the fuU EM 



qJ^^{P + Q,P)^ q^riP + Q,P) + 2^Aa2 = a^G^'iP + Q) - G^'(P)^73, 



(60) 



where Eq. (44 1 as been applied. Therefore the full EM vertex indeed obeys the important GWI in the Nambu space. 



Having the expression of full EM vertex, one can show that it also respects the g-limit GWI 

d±{P) 



hm f"(P + g, P)U=o = \ ^ ^ = <T3 



which is a sufficient and necessary condition for the compressibility sum rule [B]. 

dn 



9/1 



= -X""(w = 0,q^ 0). 



(61) 



(62) 



This is proven as the following: 



dn 

dfj, 



Tr 



V ( dGjP) 



(^3 



= -Tr 



p 

TtJ2 {t''{P,P)G{P)l°{P,P)GiP) 



^3 



9S(P) 



)G(P)a3 
0,q- 



0). 



(63) 



The g-limit GWI (61) has a profound physical implication for interacting Fermi gases. The LHS of it is associated 
with the equation of state derivable from one-particle correlation functions, while the RHS involves the response 
function evaluated by two-particle correlation functions. Hence it builds a bridge connecting the one-particle and 
two-particle formalisms. 

There is a subtlety that needs to be addressed here. It is on whether the g-limit GWI (61) can be derived from 
the GWI (551, or equivalently, whether the g-limit GWI serves as an independent constraint. Here we argue that the 
GWI imposes no constraint on the form of the g-limit GWI. The reason is because in the limit cj = and q — > 0, the 
GWI (55) becomes q • r(P, P) = limq_j.o[cr3G~'^(P + Q) — G~^(P)cr3] j^^oi where T is the spatial component of the 
vertex function. However, the g-limit GWI is an identity regarding so the GWI does not reveal any information 
about the g-limit GWI in the limit w ~ and q — > 0. Thus one should treat the g-limit GWI as an independent 
constraint of a linear response theory. The problem of how to obtain a consistent expression for the compressibility as 
discussed in Ref. [3] thus can be rephrased as how a linear response theory can satisfy the compressibility sum rule, 
or more directly whether the g-limit GWI could be satisfied. 

Now we show that the 0-th component of the F'* given by Eq.(59l satisfies the g-limit GWI. Since gi2(g) 





-g2i(g) = g§3(g) = -g§2(Q) = O when w = O, then n§(g) in this limit. Hence we have 



lim F"(P-f g,P) 

q— >0 



j=0 



(73 -ai lim H;(g) 

q— !-0 



0=0 



0-3 



Q?3(0,q 



0) 



gii(0,q^O) 



(64) 
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where 



l-2/(£;p) , df{Ep) 



Qn(0,q^O) = A^E^[ 



1 rl-2/(i?p) , df{E^) 



E^ 



dE^ 



(65) 



Note the self-energy operator is given by E = — (TiA, hence we have 

9A 



9S] 
9/i 



(66) 



We need to evaluate ^ from the gap equation (31 1. Differentiating both sides with respect to ji one gets 



Ep #2 ( 



2/(iJp 



) 8/(gp 



-Q°3(o.q^o)_i,,n;(Q) 

Qii(0,q^O) q^o 



Plug this into the right-hand side of Eq. ( |64[ ) , we prove the Q- limit GWI ( 63 1 for the full EM vertex. 
Here we summarize our key results in the Nambu space: 

• Response kernel: 



i^^-(Q) = Tr ^ (f'^(P + g, P)G{P + Q)r{P, P + Q)G{P)) + 

• Full EM interacting vertex: 

f^(P + Q, P) = r{P + Q,P)- ain5'(Q) - a2U^{Q), 

• Generalized Ward identity: 

a3G-\P + Q)-G-\P)a3^q^r''{P + Q,P), 

• Q-limit generalized Ward identity: 



t'iP,P)^a3- 



9S(P) 
5/x 



/-sum rule: 



r+oo ^2 

J-oo m 



(67) 



(68) 



(69) 



(70) 



(71) 



(72) 



We emphasize that a consistent linear response theory for BCS superfluids should satisfy the GWI, Q-limit GWI, and 
/-sum rule. 



2. WI and Q-limit WI in the one- dimensional space 

We have shown the details of the CFOP theory for BCS superfluids in the Nambu space. Instead of formulating 
BCS theory in the matrix form in the Nambu space, one may use the Green's functions (29) to formulate BCS theory. 



The Leggett-BCS theory of BCS-BEC crossover [20l |2T] follows this path and we call this "BCS theory in the one- 
dimensional space" . This representation is convenient when generalized to the theories for BCS-BEC crossover. One 
may build the CFOP theory step by step in the one-dimensional space similar to what we have done so far. Here we 
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simply extract the key results from those in the Nambu space since the underlying physics is the same. We only need 
to find the full EM vertex in the one-dimensional space such that the EM response kernel can be evaluated by 



i^^"(Q) = 2y r'^(P + Q, P)G{P + Q)Y{P, P + Q)G{P) + (73) 

^ — ' TO 

P 

where is the bare EM vertex given in Section [TTTl Moreover, the full EM vertex must obey the WI 

q^J^^{P + Q,P)^G-\P + Q)-G-HP) (74) 

and the Q-limit WI 

.imr"(P + 0,P,Uo=^^l-'5^. (75) 

The expressions of the gauge-invariant response kernels should not depend on the space in which we evaluate the 
response functions. Hence K'^'^{Q) given by Eq. (731 must be the same as that in Eq. (56). We can derive F'^ by 
using this relation. We emphasize that the two sides of the Wis (74) have very different meanings. The left-hand side, 
which involves the interaction between fermions and the external field, is a single-particle process. The right-hand 
side, however, contains many particle effects since the Green's function contains the self energy that represents the 
interactions among particles. 

The bare EM vertex in the Nambu space can be expressed as j^{P+Q, P) — diag(7^(P-|-(5, P), —j^{—P, —P — Q)). 
We define 

n^(Q) = _ni'(Q) + *n^(g), n^(Q) = ^n^iQ) - m^(g), (76) 

which satisfy 

q^n^iQ) = 2A, g^n'^(g) = -2A. (77) 



In fact these two equations are the off-diagonal terms of the GWI ( |55[ ) in the Nambu space. Hence from Eq.(|59[) the 
matrix form of the full EM vertex is given by 

F(P + g,P)-l ^^^^^ -^^{-P,-P-Q) )■ 



This type of expression was earlier obtained in Ref . [12] . Substituting Eqs.([78| and ( |29[ ) into the expression ([56|, we 
have 



K^-(Q) = -h^- + y UA'(p + p^GiP + Qh'iP, P + Q)G{P) + H^(Q)P(P + Q)Y{P, P + Q)G(P) 
p 

- r{P + Q, P)F{P + Q)r{-P - Q, -P)F{P) + H^(g)G(-P - Q)r{-P - Q, -P)F{P) 

+ n^iQ)GiP + Q)YiP, P + Q)F{P) - 7^(-^, -P - Q)F{P + g)7"(^, P + Q)F{P) 

+ W{Q)F{P + g)7"(-^ - g, -P)G(-P) + 7^(-P, -P - Q)G{-P - Qh^i-P - Q, -P)G(-P)](79) 

Changing variables by — P P + Q for the terms containing 7^(— P — Q, — P) and using F{P) = P(— P), we get 

K^'^iQ) ^ b'^iP + Q> P)GiP + Q)l"{P, P + Q)G(P) + H^(g)P(P + g)7"(-P, P + Q)G{P) 

p 

+ n'^(g)G(p -I- g)7"(^', p + Q)f{p) - j^i-p, -p ~ q)f{p + q)Y{p, p + Q)f{p)\ + -h^^r (so) 



Substituting the expression F{P) = AGo(— P)G(P) into Eq. (80 1 and comparing with Eq. ([73]), one can find the full 
EM vertex 

F'^(p-fg,p) - 7'^(p + g,p) + AH'^(g)Go(-P-g) + An^(g)Go(-P) 

-A2Go(-P)7^(-i^, -P - Q)G„{-P ~ Q). (81) 

The second and third terms can be shown to contain collective-mode effects The fourth term corresponds to the 
Maki-Thompson diagram [7]. One can verify that this full interacting vertex obeys the WI (74 1 in the one-dimensional 
space. Contracting both sides of Eq.(81 ) with g^, we have 

g,,F^(P + Q,P) = Go i(P + Q)~ G^\P) - 2S(P + Q) + 2S(P) - M^^^M^ (G^ ^(-P) - G^^-P - Q)) 



14 



= Go\P + Q)-G^\P)~^{P + Q) + ^{P) = G-\P + Q)~G-\P), (82) 

where we have used the fact = — A^G'o(— -P) for BCS superfluids. Moreover, the /-sum rule can be shown to be 

satisfied. 

From the expression ( |81[ ) we can find that the many-particle effects are indeed included consistently in the correction 
of EM vertex. The collective modes, which correspond to the poles in the response functions, are many-particle 
effects. There is a gapless mode associated with the Nambu-Goldstone mode due to the spontaneous breaking of 
the U(1)em symmetry [8]. Because the contributions of the Nambu-Goldstone modes are properly included, gauge 
invariance of our theory is restored. 

Now we verify that the 0-th component of the full EM vertex satisfies the Q-limit WI (75), from which the 
compressibility sum rule can be derived [B]. By using G^^(P) = Gq^{P) — S(P), this can be shown straightforwardly. 

= -2^r°(F,P)G(P)7"(P,F)G(P) = -/^°°(cj = 0,q^O). (83) 
p 

For non- interacting Fermi gases, we have 'y°{P + Q, P) = 1 and E = 0, hence the Q-limit WI is automatically satisfied. 
Note that limq_>o Go(P + Q)\ui=o = Go{P), we evaluate r°(P + Q, P) in the limit w = and q 

limr"(P + Q,P)Uo = l + Ahm (n"(Q)+n^(Q))|„^oGo(-P)-A2G2(-P), 

q— q— ^0 

= l-2Alimn;(Q)UoGo(-P)-A'G2(-P). (84) 

q— J^O 



Using E(P) = -A^Go(-P), the RHS of Eq. (751 is 



= 1 + 2A^Go(-P) - A'Gli-P), (85) 



where the identity 9^Go(-P) = -Gl{-P)df,GQ^{-P) = -Gl{-P) has been applied. Comparing Eqs. ^ and ([85]), 
we found that the Q-limit Ward identity holds for BCS theory only when 

o A 

— =- limnO(Q)U^o. (86) 



This has been shown in Eq. ( 67 1 so the Q-limit WI is also respected, which then guarantees the compressibility sum 
rule. 

To compare with the results in the Nambu space, we also list our key results of the CFOP theory in the one- 
dimensional space 

• Response kernel: 

K^"'{Q) = 2 y r^(P + g, P)G(P + Q)Y{P, P + Q)G{P) + -h^'\ (87) 



• Full EM interacting vertex: 

r^(P + g,P) - 7^(P + Q,P) + An^(Q)Go(-P-0) + An^(Q)Go(-P) 
-A2Go(-P)7^(-P, -P - Q)Go{-P - Q), 

• Ward identity: 

q^.r^{P + Q,P) = G-\P + Q)- G-\P). (89) 

• Q-limit Ward identity: 

r"(P,P) = i-^, (90) 
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/-sum rule: 



+00 



dujuJXpp{uj,q) = n- 



(91) 



The WI (or GWI) of the CFOP linear response theory for BCS superfluids guarantees that it is gauge invariant. 
In Appendix [E] the explicit gauge invariance of the CFOP theory is studied from another point of view based on a 
"generalized gauge transformation" . 



C. Review of Nambu's Linear Response Theory 



Since we have the full EM vertex f from the CFOP theory, it would be helpful to compare it with the results 
from Nambu's integral-equation approach [i]. A parallel discussion for relativistic BCS superfluids can be found in 
Ref. fl4'. In BCS theory of conventional superconductors, the self energy is approximated by an integral equation 
which consists of a ladder approximation for the electron-phonon interaction. Nambu proposed that the EM vertex 
should be corrected in the same way as that for the self energy. Hence the EM interacting vertex follows an integral 
equation 



f ^(P + Q,P)= T{P + Q-P)+9Y. '^3G{K)f^'{K + Q, K)G{K + Q)a^. 



(92) 



K 



If a vertex is a solution to this equation, it automatically satisfies the GWI ( 55 ) . We give our own proof in Appendix [C| 
Ideally, if we know how to solve this integral equation, we can further calculate the gauge-invariant response kernel 
K^^ by Eq. (56). Unfortunately, very little is known about the solution. In general, one may expand the RHS of 
the equation as a series of g by the iteration method, and then truncate it at some order. However, this will not 



produce a gauge-invariant solution. Moreover, the integral equation (92) is a vector equation while the GWI (55) is a 



scalar equation so they have different degrees of freedom. This suggests that there should not be a rigorous one-to-one 
correspondence between the solutions to the integral equation and the EM vertex respecting the GWI. As pointed 
out by Nambu [1], the integral equation is not only consistent with the GWI associated with the EM vertex but also 
consist with the GWIs associated with three other interaction vertices or gauge transformations (as shown in Eq.(4.4) 
of Ref. [4 ). 



In fact, what is equivalent to the GWI (55) is the contracted integral equation given by 



g^f + Q,P)^ q^riP + Q, P) + gY.^zG{K)q^V^{K + g, K)G{K + Q)a3. 



(93) 



K 



The sufficient condition of this proposition is that the contracted integral equation (93) can lead to the GWI, which 



is proven in Appendix [Cj The necessary condition of this proposition is that any EM vertex obeying GWI must also 
satisfy Eq. ( 93 ) , but not necessarily the integral equation ( 92 ) . Importantly, this is pointed out by Nambu in his 
seminal paper [4 and by Schrieffer [THj. We briefly outline the proof here. Substituting Eq. (55) into the RHS of 



Eq. ( 93 ) , we have 



RHS = g^r(^ + Q,n+.9l]'^3G(i^)(^T3G-\X + 0)-G-\X)a3)G(if + Q)a3 



K 



= q^.r{P + Q,P)+ 2iAa2 = q^T'^iK + Q,K)= LHS, 



(94) 



where Eq. (CI ) has been applied. Therefore, any gauge invariant EM vertex F'^ (including the vertex from the CFOP 



theory) must be a solution of the contracted integral equation ( 93 ) , but not neces sari ly a solution of the integral 
equation (92). Furthermore, F''' may differ from the vertex F'^ obtained from Eq. (92) by a gauge transformation 
p'M = _|_ T^M^ where satisfies the Lorentz equation q^x^ = 0. Such a gauge transformation may not necessarily 
be expressed as x'' — J , where / is a matrix in the Nambu space whose elements are harmonic functions. Another 
example is given by = Y['^{Q)C^ where C is an arbitrary constant matrix in the Nambu space with at least one 
nonzero element. H^((5) is given by Eqs. (57) and q^IiilQ) — 0, hence j<f satisfies the Lorentz equation. Under the 
gauge transformation F'' — ^ F^ 
CFOP theory transform as 



X^ , the 0-th component of the vertex and the density response function from the 



OA 1 

fO(P,P) ^ f"(P,P) - —C, Xpp(Q) ^ Xpp(Q) - -ImTr^ {lVl{Q)CG{P + Q)a^G{P)) . 



(95) 
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Then one can verify that the Q-hmit GWI (71) and the /-sum rule (72 1 can not be satisfied simultaneously under 
the above transformation even though the GWI ( 70 1 is respected anyway. 

Here we have two remarks on the consistency of a linear response theory for BCS superfluids. Firstly, from the 
proof given in Appendix [C] Nambu's approach tells us that if the self-energy of an interacting Fermi gas satisfies 
Eq.(Cl), then the vertex given by the integral equation (92) must satisfies the GWI. However, the theory of BCS 
superfluids is not the only theory that satisfies Eq. (CI ) and as shown in Ref. [4j there are theories with other types 
of symmetries that satisfy Eq. (CI ). Hence Nambu's vertex given by Eq. ([92]) may not be specific to BCS superfiuids 
and it can be more general. Secondly, the Q-limit GWI and /-sum rule are not directly linked to the GWI. They 
should be imposed as separate constraints for a consistent linear response theory. The case of non-interacting Fermi 
gases is special because the bare EM vertex ^^{P + P) is already the full vertex and no correction is needed. If one 
formulates a linear response theory and finds a gauge-invariant vertex which obeys the GWI, one can further calculate 
the response kernel by Eq. (68 1. However, this response kernel may not satisfy the compressibility sum rule and /-sum 
rule. Therefore, as we emphasized, the GWI, Q-limit GWI, and /-sum rule should be independent constraints for 
a consistent linear response theory for BCS superfiuids. We cannot fully check whether Nambu's integral-equation 
approach satisfies all these criteria because finding an exact solution to the integral equation is the bottleneck, but 
the CFOP approach does satisfy all those constraints as we have demonstrated. 



D. Comparisons between different linear response theories for BCS superfluids and Meissner Effect 



Now we compare the two vertices given by the CFOP theory and Nambu's integral-equation method. In Ref.[l], 
Nambu managed to solve Eq. (92 1 with the aid of Eq. (551 under the following conditions: (1) only the zeroth order 
of g is considered, (2) w and q are both small, and (3) it is at zero temperature. Here we compare the result ( |59[ ) with 
Nambu's under the same condition. Taking the same limits, the response functions Q|i, Q12 and Q^a vanish because 



of the particle-hole symmetry. Therefore n°((3) = — H*'((5) 



2{Q) 



and n(Q) = -n(Q) 



expressions of the response functions, one can verify that (see Appendix |D[) 



■ Q23(Q) 
'Q22(Q)' 



.N{0)l 



Q23(Q) 



N{0)cl 

I \ q, 



2A ' A 
where iV(0) is the density of states at the Fermi energy and c 



2A2 



From the 



(96) 



1 fcf 

^3 m 



-^vp is the speed of sound. Here kp is the 



k 

Fermi momentum defined by n = Therefore we have 



and 



n°(g) = -nO(Q) 



r"(P + Q,P) = (73 + 21(72 



Ao 



n(Q) = -n(g) 



2Ac? 



ciq^ 



Auj 



T{P + Q,P) 



2i(72 



2 ~ c^q'^ 



(97) 



(98) 



cj — Csq is the dispersion of the gapless collective mode. These results are exactly the same as those found by Nambu 
[4] by taking the same limits. 

Since the CFOP vertex and Nambu's vertex both satisfies the GWI, they are at most off by a term with (/^x^ = 0. 
Xn can be expressed by a harmonic matrix / in the Nambu space of the form = 9^/. In momentum space, it is 
Xii = ^'Zm/i which vanishes as q^ — !■ 0. Therefore, at zero temperature and in the low frequency and momentum limit, 
the CFOP linear response theory agrees with Nambu's approach but in general they can be different. 

Before closing our discussions on the density channel, we remark that the collective modes do not contribute to the 
Meissner effect, where one can show that a finite superfluid density leads to perfect diamagnetism [5]. This remark 
justifies the standard calculation of the Meissner effect and superfluid density [2], where one ignores the fluctuations 
from the order parameter and still obtains the correct results. Although one may use a fully gauge- invariant linear 
response theory to demonstrate this result |12j . here we directly evaluate the collective-mode contribution in the 
Meissner effect. Following Ref. [5], the Meissner effect is associated with the current-current correlation functions, 
which can be inferred from the transverse components of the response kernel. Here we briefly sketch why the collective- 
mode effects do not contribute to the transverse components of K'^^{{)^q) as q — > 0. From Eq. (41), the response 
kernel K^^ is given by 



Q11Q32Q23 



Q22Q3lQi3 — 2(5i2Q3iQ2 



Q11Q22 — Q12Q2I 



(99) 
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where the second term is associated with the collective modes. A tensor P*^ can be decomposed into the longitudinal 

•<-> <->■.. 

and the transverse parts Pl and Pt, where Pl = q • P • q, Pt = Q2iP" — Pl)/'^, and q is the unit vector along 
q. Assuming that q is parallel to the z-axis, in the limit q — )■ one can show that and Q32 start with linear 
dependence on q. Therefore limq_^o Qa; • Qaj = limq_j.o q • QsiQsj • q- When one evaluates the transverse components 



of Eq. (99), all of the collective-mode terms cancel in the limit q — > so the demonstration of the Meissner effect is 



insensitive to the collective modes. 

V. SPIN LINEAR RESPONSE THEORY OF BCS SUPERFLUIDS 

We now formulate a generalized spin linear response theory which is similar to its counterpart in the density channel. 
Using the notation of the Nambu spinor (251, the Hamiltonian ^ with the BCS approximation in momentum space 
is given by 



pq 



Here the fluctuation of the order parameter Aq is in the spin channel and should not be related to the fluctuation 
in the density channel. We follow the same procedure as what we did in the density channel. The order parameter 
is separated into two parts with one denoting its equilibrium value and the other denoting the perturbative part. 

Furthermore, by introducing the spin interacting vertex 7g (P + Q,P) = 7s (p + q, p) = (1, ^'^^'^3) in the Nambu 
space, the Hamiltonian can be expressed as H = Hq + Hg, where 

i/0 = E*P^P*P' ^S=E*I.+q(^lq^l+^2q^2+A^q7s^(p + q,p))vI/p, (101) 

p pq 

Similar to its density counterpart, the interacting Hamiltonian can also be cast into the form 

^s = E*l+q'^q ■^s(p + q,p)*p (102) 
pq 

by introducing the generalized driving potential and generalized interacting vertex 

*q= (Alq,A2q,A^q)^, Eg (p + q, p) - (^1 , ^2 , 7s (P + P) ) ^ , (103) 

The Heisenberg operator is defined as 0{t) = e^'^Oe^^'^ . Hence, when the external field is weak, the generalized 
perturbed current in the spin linear response theory is given by 



<5 Js(t, q) = V(vl/t (r)Ss(p + q, P)vl'p+q(r)) + -S^^h^^^ A^ir, q). (104) 



Here SJg^ denotes the perturbed spin current and SJ^-^ ^ denote the perturbations of the gap function. The spin linear 
response theory can also be written in a matrix form 

(5Js(w,q) = Qs('^,q) ■ *('^,q) 

Qsii(w,q) Qsi2('^,q) (5gi3(w,q) \ / Ai(a;,q)\ 

Qs2i(w,q) Qs22(^,q) Qs23('^.q) A2(a.,q) . (105) 

Qg3i(w,q) Og32(^,q) Q^;,iu;,q) + ^h^- J \A,iu;,<Dj 

Following the same steps as what we did previously, the spin response functions are also expressed by 

Qs^J (c^, q) - TrT ^ ^ (Eg, (P + Q, P)G(P + Q)Es, (P, P + Q)G(P)) , (106) 

icjn p 

The expressions of these response kernels are given in Appendix |Aj One can verify that Qsii = Qii, Qsi2 ~ Q12, 
Qs2i = Q21 and Qs22 = Q22 so the block of the matrix Qsij for the fluctuations of the order parameter is exactly 
the same as that in the density channel. It can also be shown that Qg^^ = Qs23 ~ Qssi ~ Qs32 ~ 0. Hence 
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the perturbation of the order parameter decouples from the perturbation of the spin current. By applying the self- 
consistent condition SJsi,2 = —^^1.2, we obtain 

Qsii(w,q) (5si2(t^,q) \ /Ai(a;,q)\ 
= I Qs2i(c^,q) Qs22(^,q) A2(aj,q) , (107) 

Qr(^.q)/ \M^,ci)J 

where Qsii = |+Qsiii Qs22 = | + '9s22 and Qg'' — Qg^^ + ^h'^'^ . The perturbations Ai and A2 in the spin channel 
can be further shown to be zero, which are very different from their counterparts in the density channel. The spin 
response function becomes 

SJ^^iu,q) = (w,q)A,(c^,q). (108) 

We emphasize that the fluctuations of the order parameter automatically decouple from the spin response and there 
is no contribution from the collective modes in the spin response function. The mechanism behind this decoupling is 
that the BCS order parameter does not break the 11(1)^ symmetry. Even though one treats the spin linear response 
in the same way as what we did for the EM response, the unbroken U(l)2 symmetry leads to a significantly different 
result. 

The U(l)z invariance of the linear response theory is satisfied by the GWI 

gMQr(Q) = 0, (109) 

which leads to the conservation of the spin current q^SJ^ — 0. Before proving this statement, it is important to notice 
that the spin interacting vertex obeys the GWI associated with the V{l)z symmetry 

q,%{P + Q,P) = G^\P + Q)- Go'iP) = G-\P + 0) - G-\P). (110) 



This leads to the GWI ( 109 1. The second equation is due to the fact that the BCS self energy in the Nambu space is 
give by S = — Acti and is independent of the four-momentum. The proof of the conservation of the spin current is 
similar to the derivation shown in Eq. ( 47 1 . 

The spin susceptibility can be evaluated from the spin response kernel by xss = ^"I^iQggg. The /-sum rule can 
be shown to be satisfied in the spin channel. 

+00 ^2 

da;a;xss(w,q) = n — . (Ill) 
-00 m 

Following the same argument in the density channel, we have 

- - / duujlniQll^ ^ dulxn(LoQl%) = dc.Im(q • Qggg). 

^ J —00 ^ J —00 ^ J —00 



Comparing Eqs.(A19l to (A13), we see that Qg33 — Q3I. Hence by using the lemma (50) we get 

- - / dwwImQggg = / (ia;Im(q • Ql^a) = — . (112) 

J —00 J —00 ^ 

In the one-dimensional space we follow the same steps in the density channel to find the U(l)z gauge-invariant spin 
vertex Tg^ which satisfies the WI ([74|. The spin response kernel can be expressed by 

n 
m 



(Q) -EE ^sAp + p)G{p + Qhup^ p + Q)G{p) + -h^\ (113) 



where the bare spin interacting vertex is given in Sec. |III| Im portantly, Eq.(113l should give the same expression 
for the spin response function as the 33-component of Eq. (llOGf when the h^^ term is included 



Qr(Q) = Tr E (7s + P)G{P + QWsiP. P + Q)G{P)) + "^h^r (114) 



p 

Note that the spin interacting vertex in the Nambu space can be expressed as 



%iP ^Q.P)-{ ''^^^'^ t ,3^,(-P, - Q) ) ■ (11^) 
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Substituting Eq. (1151 and the propagator (28) into the expression (114), we get 

= ^h^- + Tr^ (7g^^(P + Q, P)G{P + Q)Ys^P, P + Q)G{P) + ^>^^{P + Q, P)F{P + g)7si(-^ - Q, -P)F{P) 
p 

+ 7s^(-^, -P - Q)F{P + Qhs^P. P + Q)F{P) + T^^C-^', -P - Q)G{-P - QHii-P - -p)G{-p)) 

= ^h^'- + TtJ2 {istiP + P)G{P + QhkiP^ P + Q)G(P) + 7s^(-^, -P - Q)^(P)7si(^, P + Q)F{P + Q) 
p 

+ 7s^;(-^^, -P - Q)F{P + 0)7st(^, P + Q)F{P) + 7^^(^ + Q, P)G{PhUP, P + Q)G(P + Q)), (116) 

where we have changed variables by — P — > P + Q for the terms containing 7g^(— P — Qi^P) and applied P(P) = 
P(— P). Using P(P) = AG'o(— P)G(P) again and Eq.(113l, one can see that the following expression guarantees that 
the expressions for Q'^'^iQ) are the same. 

rg,(P + Q, P) = 7s^(^ + Q,P) + A^Goi-PHA-P^ -P ~ Q)Goi-P - Q) (117) 

for (T |. The second term corresponds to the Maki-Thompson diagram [7), which is the same as its EM counterpart 
except there is a minus sign. Now we show that this full spin interacting vertex satisfies the Ward identity (74). Using 
Scr = —Sg- and Eq.([l5]), we have 

qA(P + Q^P) = Sa{Go\P + Q)- G^\P)) - S,A^Go{-P)Go{-P - Q){G^\-P) - G^\-P - Q)) 

= S^{G^\P + Q) - G^'iP)) - S,{j:{P + Q) -^{P)) ^ S4G-\P + Q) - G-\P)). (118) 

Finally we summarize the central results in this section. In the Nambu space, we have 

• Response kernel: 

Qr(Q) = Tr ^ (7^(P + Q, P)G(P + Q)7s^(P, P + Q)G(P)) + (119) 
p m 

• Spin interacting vertex: 

7s^(P + Q,P) = (l,Pila3), (120) 



• Generalized Ward identity: 

q^l^iP + Q,P)^ G-\P + g) - G-\P). (121) 

In the one-dimensional space, we have 

• Response kernel: 

(Q) = E E ^s.{P + Q-^ P)G{P + Q)iUP: P + Q)G(P) + (122) 

• Full spin interacting vertex: 

rg^(p + g, p) = i^,{p + Q,P) + A^Goi^PH^i-p, -p - Q)G„{-p - g), (123) 

• Ward identity: 

Q,^^sAP + Q,P)^ S,{G-\P + g) - G-\P)). (124) 



In both spaces, the /-sum rules are satisfied 



+00 ^2 

c^wa;xss('^,q) = " — ■ (125) 
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Below Tc, the U(1)em symmetry is broken by the condensed Cooper pairs while the U(l)z symmetry remains intact. 
Hence the two linear response theories produce significantly different results below Tc. Importantly, the collective 
modes coming from the breaking U(1)em symmetry only couples to the density response function but decouples from 
the spin response function. This indicates that the difference between the density and spin channels arises only in the 
presence of a condensate that only breaks the C/(1)em symmetry. Above Tc, both U(l) symmetries are respected since 
there is no Cooper-pair condensation. The effective Lagrangian after the BCS approximation is identical to that of a 
non- interacting Fermi gas when the gap A vanishes. By dropping Qn and Q2i which are not defined above T^, one can 
verify that the response functions of the two linear response theories give the same result above Tc. In the presence of 
pairing fluctuation effects, the amplitude of the pairs are not necessarily the order parameter since finite-momentum 
pairs may coexist with the condensate of Cooper pairs [21] . The difference between the two response functions could 
be shown to be still valid and one can use the difference between the density and spin structure factors to construct 
a quantity similar to an order parameter for detecting the phase coherence of atomic Fermi gases [5] . 



We have shown that the CFOP approach of the linear response of BCS superfiuid is a computational manageable 
scheme that satisfies important constraints including Ward identity, / sum rule, Q-limit Ward identity, and compress- 
ibility sum rule that guarantee charge conservation and a consistent expression for the compressibility. The CFOP 
formalism provides a paradigm for studying linear response theories in interacting many-body systems in the presence 
of spontaneous symmetry breaking. The spin linear response theory complements the story of the CFOP theory and 
demonstrates the different roles played by the collective modes in the superfiuid phase. Going beyond mean-field BCS 
theory requires considerations of non-condensed pairs and there have been different approaches [7j . We emphasize that 
linear response theories of those beyond-BCS theories should be subject to the same consistency constraints discussed 
here. In addition to conventional superconductors |12j . our formalism may be useful in the study of ultra-cold atoms 
[SJ m] and nuclear physics [531 HI] where linear response theories of BCS superfluids are frequently implemented. 
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Appendix A: Detailed expressions for response functions 



The following are the response functions for non-interacting Fermi gases: 




(Al) 




(A2) 




(A3) 



The following are the EM response functions of BCS superfiuids from the CFOP theory: 



Qn(c^,q) = I] [(1 + 



EpEp 



A2 



) 




•^-{E++Ep) 



-[l~f{E+)~f{E^)] 



p 




-A' 



2 



Et -E, 



(A4) 
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<3i2(w,q) = -(921 (w,q) = -i^^ [ 



^+ £^ i-f{E+)-fiE^) g+ /(ig+) -/(£;-) 

'^^^E+ Ep' io^-{E++Ep)^ ^E+ Ep' uj-^-{E+ -Epf 



, (A5) 



Q;3(c.,q) = g°i(a;,q) = A^ 



e+ + liE+ + E-)[l - f(E+) - f{E-)] {E+ - E-)[f{E+) - f{E-)] 



Ep Ep 



cj^-{E++Epy 



+ 



uj^-{E+-Epy 



, (A6) 



Qi „^ _ OM.; V ^'^ ^ {E^ - E-)[l - fjE^) - fjE^)] + )[/(£+) - KE^)] 

^13\'^'^) — ^31\^^W — 2^ — — L — s :;3T7; r 



mE+Ep 



u^-{E+ + Ep)^ 



:2-{E+-EpY 



(A7) 



<922(w,q) 



-[f{E;)-f{E- 



E+Ep 



P-{E+-Epr' 



(A8) 



E+Ep 



{E+ + i?-)[l - f{E+) f{E-)\ ^ {E+ - E^)[f(E+) f{E^ 



u^-{E++Ep)^ 



LO^-{E+-Epy 



, (A9) 



Q23('^>q) = -Q32('^>q) = 



p^C+-^p r(i^+ + i?-)[i -/(£;+) ^ )[/(£;+) -/(i^p 



m E+E, 



p ^p 



;2 )2 



+ 



E+Ep ^w^-{E+ + Ep)^ 



E+ + E- 
P P [1 - f{E+) - f{E-, 



£■+ - E- 
" :[f{E+)~f{E-, 



^-{Et-Ep^ 



(All) 



-(1 + 



Jt^p+^. 
EpEp 

5p2P__Z_'\ 
E+Ep 'c 



:^-{Et+E-Y 



[1-/(E+) -/(£-)] 



K-Ep 



Ep)' 



[f{E+)-f{E-)] 



(A12) 



Q^3('^,q) = Q3°3('^,q)='^E^[( 



Ep^ uj^-{E++Ep)^ ^E+ Ep ^ uj^ - {E+ - Ep)^ 



(A13) 



The following are spin response functions of BCS Superfluids following the same structure of the CFOP theory 



—d _ ^p^p \ p p 



i;+£;p ^a;2_(£;+_^-)2 



[f{E+) - f{E-^ 



(A14) 



<3si2(w,q) = -(5s2i(w,q) = -i^E [ 



e-,l-fiE+)-f{E-) ,e Cn/«)-/(^p) 



Ep^ uj^-{E+ + Ep)^ ^E+ Ep^uj'-{E+ -EpY 



, (A15) 
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Cp^p+^\ 5t_Ep [/(^+)_/(^^)]\ (A18) 

E+Ep ' uj-^ - {E+ - EpY^^ ^ ■'^ P^^i' ^ ' 



^S33l -4S331 Z^^\V^+ i;p^ ,^2_(^++^p)2 V ^+ ^ ^2 _ (^+ _ )2 / ' ^ ^ 



P 



+ "n^ K^ -%i ~ 

Here we outline the proof of Eq. ( |A17[ ). When the index = 0, from Eq. (106) we have 



Q0i3(c.,q) = Y,{G{P)F{P + Q)-G{-P-Q)F{P) + G{P + Q)F{P)-G{-P)F{P)). (A21) 
p 

Changing variables by P — > — P — Q and using the fact F{P) — F{~P), the second and fourth terms inside the 
bracket become G{P)F{P + Q) and G{P + Q)F{P) respectively, which cancel the first and third term inside the 
bracket respectively. Hence we have Qgi^{uj, q) = 0. Similarly, if the index fj, = i, we have 

Qsi3(w,q) = 5]^^(G(P)F(P + Q) + G(-P-Q)P(P) + G(P + g)F(P) + G(-P)F(P)). (A22) 
p 

Here we change variables by P — > — P — Q again so the pre- factor - — - — Hence all terms inside the 

bracket cancel out and we conclude that Qsi3(t^, q) = 0. Following the same steps, one can prove that Qg23(i^, q) = 
<3s3i('^'q) = Qs32('^'q) = 0- 

Appendix B: Proof of the Lemma ( |50[ ) 

From the expression shown in Appendix |Xj we have 



Ep Ep Ep + Ep V cj — Ep — Ep ui + Ep + Ep 

Ep' E+-Ep \uj-E+ + Ep uj + E+-EpJi 



Hence 



+ 00 ^ 

da;-Im[q- Q°3(a;,q)] 
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i-oo <• -Ej E^' E*+Ef V> p pi V P pV 



V 77'+ ip— / rp-r- 771— P ^ P P' 



<5(w - + ) - S{u; + E;- iJp )) } 



^ ^ [1 - - ))] - 2E ^ [1 - ^(1 - 2/(i^p+))] . 



(B2) 



We change variables by p ^ p + ^ in the first term, and change variables by p — > p — ^ in the second term to get 



+ CXD 



da;-Im[q • Q^giw, q)J = ^ ^ ^ 



l-|^(l-2/(i?p)) 



(B3) 



where the number equation (31 ) has been used. This proves the lemma. 



Appendix C: Integral Equation of EM Vertex and GWI 



Here we prove that the vertex determined by the integral equation ( 92 ) obeys GWI ( 55 ) . We will use the following 
equality 



g J2 (^3G(P) - G{P + Q)<T3) = -(^sS - Sag). 
p 



(CI) 



To prove the proposition, we only need to show 

asG-^P + 0) - G-\P)a3 - q^riP + Q, P) + 9Y,cj^G{K)q^t^^{K + Q, K)G{K + Q)a^, (C2) 



K 



which is equivalent to 



(fTgE - Sda) = g ^ a:iG{K)q^J:^'(K + Q, K)G{K + Q)fT3. 



K 



From E = Gg - G~^{P) one concludes that 



G{P)G^\P) = 1 + G{P)±, Gq i(P)G(P) = 1 + SG(F). 



(C3) 



(C4) 



Now we turn to the prool of Eq.(C3). By substituting Eq.(55) into its RHS and repeating the process, we get an 
iterative equation 



RHS of Eq. (C3) 



gJ2^3G{P)q^riQ)G{P + QW3 
p 

V '^3G{Pi)a3G{P2)q^r{P + g, P)G{P2 + Q)<J3G{Pi + Q)a3 + ■■■ 
P1P2 

00 

Y^g' '^3G{Pi)---a3G{P,)q^j^'{P, + Q,P,)G{P, + Q)a3---G{Pi+Q)a3 

i=i Pi--Pi 

00 i i 

Y l[['^^G{Pk)]qf,J'^{P^ + Q,P^)l[[G{P,+^-k + Q)'J:i\, 

i=i Pf-Pik=i 



(C5) 



After inserting the Ward identity (54) for the bare EM vertex and using Eqs.(C4), we get 



RHS of Eq.(|C3| 

00 

Yd" Y ^3G{Pi)---o^G{Pi)a^G^\P, + Q)G{P, + Q)a^---G{Pi+Q)a:i 
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- E^'' E '^3GiP^)■■■a3G{P,)GoHP^)a3G{P, + Q)a3■■■GiP^+Q)a3 

i=i Pi-Pi 

= gJ2i<y3G{P)-G{P + Q)<j^) 

P 

oo 

+ Y.9"Y. '^3G{Pi)---a3G{P,)G{P,-i+Q)a3---G{Pi+Q)a3 

i=2 Pi-Pi 
oo 

+ E f^3G(Pi)---a3G(P.)^3SG(P,; + g)a3---G(Fi+Q)(T3 

i=i Pi-Pi 

oo 

- ^.9' f^3G(Pi)---^T3G(P,-i)G(P, + Q)(73---G(Pi+Q)a3 

i=2 Pi •••Pi 

C30 

- J2 fT3G(Pi)---a3G(P0Sa3G(P, + g)a3---G(Fi+Q)a3. (C6) 
i=i Pi^^^Pi 

By changing the dummy index i — > i + 1 in the second and fourth summations, we get 



RHS of Eq.(C3) 



00 i i 

E n hG(Pfc)]g^ [a3G(P,+i)-G(P,+i+Q)f73] [] [G(P,+i_,. + Q)a3] 

oo i i 

E5' E n hG(Pfe)] hS - tas] n [G{P^+l^k + Q)a3] 



(crgS - ScTg) = LHS of Eq.(C3), (C7) 



where Eq.(Cl ) has been apphed. Therefore we have proved that any vertex that satisfies the integral equation must 



also satisfy the Ward identity and hence must be gauge invariant. 



Appendix D: Evaluations of the Vertex 

From the expressions given in Appendix |Xj at T = we have 

where in the fifth hne we have used — 2m(^p + fi). Note that fi c^^ — Thus 

Similarly, for the temporal component of Qt^^ we have 

dh ^ ^ 

(27r)3 2£;2 -AEl ~ * 2A 

Therefore n°(u;,q) = ^2 ^2 and the temporal component of the full vertex is given by 

fO(P + g,P)~a3 + 2ia2^^^^. (D4) 

For the spatial component, we have 

Q23(a;,q)^*Aq-y (^^^Z^--3;^q-iy (2;o^2^^-3A^r^- 



. / d3p Ao; 2E^ N{{))uj 
g23(w,q)~* / — ^— ^— ^ = — , (D3) 



25 



(D5) 

Using ^ = ^c1, we have Q23(w, q) -i^^^^q , so 



Zl\C 

n(c.,q)^^ ^q. (D6) 



The spatial component is then given by 



Appendix E: Physical interpretation of the gauge- invariant linear response theory 

In the main text we have seen that the Wis for response functions are indeed satisfied even without knowing the 
exact form of the full EM vertex. Here we pay some attention to the gauge invariance of the CFOP linear response 
theory from the point of view of a gauge transformation for the BCS Lagrangian. In real space, the Lagrangian 
density following the BCS approximation is given by 

Cbcs = ^-^^^ - i^-^^ - /^)^3 - + Adi)*, (El) 

where ^ (as, - '^). This La grangian density is obviously not invariant under the infinitesimal gauge transformation 
^I' — > {l — ia^x)"^, ^'^(l + «o'3X) a.nd — > A^ + 9^x if the fluctuations of the order parameter are not considered. 

Now we split the order parameter into its equilibrium and perturbative parts as A — ^ A + A', where the equilibrium 
value is A and the perturbation is A'. Therefore in the CFOP theory, the Lagrangian density in real space becomes 

/:cFOP = /:bcso + C' = - - ^l)a3 + Aai)^ - ¥ {Aai + Ascia + ^,7'^)^', (E2) 

Here the subscript "0" denotes the part in equilibrium. The gauge transformation of 5" and ^'^ leads to the fluctu- 
ations of the order parameter 6A1 = and 5 A2 = — 2xA. Therefore, the following generalized infinitesimal gauge 
transformation leaves the Lagrangian density (|E1[) of the CFOP theory invariant: 



* ^ (1 - icTsx)*, ^ ^-^(1 + icr3X), A^A, 

A^^A^+d^X, Ai^Ai, A2-^A2-2Ax. (E3) 
Under this generalized infinitesimal transformation the two parts of the Lagrangian density transform according to 

-Cbcso ^ /^Bcso + *^a^X7''* - ix*^A[CTi, (73]* = Cbcsq + *^a^X7^* - 2x*^Acr2*, 

£' ^ £' + 2x*^Act2* - ^'^9^X7"*- (E4) 



Therefore £bcs is indeed invariant under the generalized infinitesimal gauge transformation ( E3 1 . It is the gauge 
transformation of A2 that compensates for the effects associated with the Cooper-pair condensation and leads to the 
gauge invariance of the CFOP theory. The Noether current associated with this generalized gauge transformation can 
be deduced by introducing the "generalized gauge space" . We define the space where the generalized external potential 



and generalized interacting vertex (see Eq. (34)) live as the generalized gauge space. Explicitly, the perturbative 
Lagrangian density is rewritten as -"^/^^^ ■ S^, where • denotes the inner product in this generalized gauge space. 
The generalized gauge transformation of the generalized external potential is 

^ ( ' \ 

or in the momentum space $ + — 2Ax . (E5) 
We define the generalized external momentums q = ( 0, 2iA, g^^)^ and q = (0, 2iA, ~q^)^ in the generalized gauge 




space. Then the generalized gauge transformation ( E5 1 can be written as 

* * -I- iqx. (E6) 
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The GWI (44) can be expressed as 



a3G-\P + Q) ~ G-\P)a3 ^ 



(E7) 



In fact, this is the generahzed Ward identity associated with the generahzed gauge transformation (E5) in the gener- 
ahzed gauge space. 

Next we address the conserved current associated with this gauge transformation. Using the self-consistent condition 
(5Ji^2 = — |Ai_2, the current in Eq. (361 can be written as 



Qii Qi2 Qi3 
Q21 Q22 Q23 
Q31 Q32 Q33 

We then define the generahzed current J = (0,0, J'^)'^ and three generahzed response-function vectors 





(E8) 



Qi 





Then the current equation ( E8 ) becomes 



Q2 



J = (Qi,Q2,QD *> 




The GWIs ( 42 1 for the response functions can also be written as 

q^-Qj = 0, fori = 1,2, 3. 
Thus the GWIs directly lead to the conservation of the generalized current 

• J = (q^ • Qi , • Q2, • Q^') • * = 0. 



(E9) 



(ElO) 



(Ell) 



(E12) 



This gives the conservation law of the EM current = 0. Therefore J is indeed the Neother current associated 

with the generalized gauge transformation. Moreover, by noting that the GWIs (Ell ) in the generalized gauge space 
can be written as 



(Qi,Q2,Q^')-q-0. 

Under the generalized gauge transformation ( |E6[ ), the generalized current transforms as 

J = (Qi, Q2, Q^) • * ^ (Qi, Q2, Q^') • (* + *qx) = J. 

Thus the generalized current is indeed invariant under the generalized gauge transformation. 



(E13) 



(E14) 
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